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ABSTRACT 



Stellar collisions are an important formation channel for blue straggler stars in globular and old open clusters. Hydrodynamical 
simulations have shown that the remnants of such collisions are out of thermal equilibrium, are not strongly mixed and can rotate very 
rapidly. Detailed evolution models of collision products are needed to interpret observed blue straggler populations and to use them to 
probe the dynamical history of a star cluster. We expand on previous studies by presenting an efficient procedure to import the results 
of detailed collision simulations into a fully implicit stellar evolution code. Our code is able to evolve stellar collision products in a 
fairly robust manner and allows for a systematic study of their evolution. 

Using our code we have constructed detailed models of the collisional blue stragglers produced in the Af-body simulation of M67 
performed by Hurley et al. in 2005. We assume the collisions are head-on and thus ignore the effects of rotation in this paper. 
Our detailed models are more luminous than normal stars of the same mass and in the same stage of evolution, but cooler than 
homogeneously mixed versions of the collision products. The increased luminosity and inefficient mixing decrease the remaining 
main-sequence lifetimes of the collision products, which are much shorter than predicted by the simple prescription commonly used 
in yV-body simulations. 

Key words. Stars: formation, blue stragglers, open clusters and associations: general, M67 



1. Introduction 



Blue stragglers are stars that appear above and blueward of 
the main sequence turnoff in the colour-magnitude diagrams 
(CMDs) of star clusters, with mass es larger than t hat of a turnoff 
star. They were first identified by Sandage (j 1 953b in the globu- 
lar cluster M3 and soon aft erwards in the old open cluster M67 
dJohnson & Sandagel 119551) . Blue stragglers have since been 
found in other globular clust ers and in open clusters of all ages 
dAhumada & Lapassetl 1 1995b. Various me chanisms have been 
proposed for their formation (IStrvketll9 93). Currently, there are 
three accepted formation channels: mass transfer due to Roche- 
lobe overflow in binary systems, and stellar mergers, either due 
to dynamical collisions or through coalescence of close binaries. 



iHills & Davl Jl976) first showed that physical collisions be- 
tween main-sequence stars are likely to occur in the dense cores 
of some globular clusters. In most environments collisions be- 
tween single stars are very rare, but binary systems can sig- 
nificantly enhance the rate of stellar collisions because of their 
much larger collisional cross sections. Stellar collisions are thus 
an important formation channel for blue stragglers even in old 
open clusters such as M67, as is demo nstrated by direct A^-body 
calculations dHurlev et alJl200U 120051) . These simulations indi- 
cate that, while all the above-mentioned formation mechanisms 
operate in different regions of the cluster and all are needed to 
reproduce the observed blue straggler population, all formation 
paths - including binary mass transfer and binary coalescence - 
are strongly affected by the dynamical evolution of the cluster. 
The blue straggler population therefore contains important in- 
formation on the dynamical history of a cluster. Extracting this 



information requires understanding not only how blue stragglers 
are formed but also how they subsequently evolve. 

In .V-body models collisions between two main-sequence 
stars are usually approximated by assuming that the stars merge 
without mass loss and replacing the end product by a nor- 
mal (evolved) main-sequence star with its ag e reset in accor- 
dance with th e assumption of complete mixing (fTout et alJ 19971 
iHurlev et al.l 120021) . Smoothe d particle hydrodynamics (SPH) 
simulations of stellar mergers dLombardi et al.ll9 95, 1996) show 
that some mass is lost during the collision and that the collision 
products are generally far from being fully mixed. To understand 
how this affects their further evolution and the predicted blue 
straggler population, detailed stellar evolution models of the col- 
lision products are needed. 

In early attempts to model the evolution of stellar merger 
remnants simplifying assumptio ns were made: the two star s 
were assumed to mix completely (Bailyn & Pin sonneaultll995h . 
a specified ch emical profile was imposed on the remnant 
(ISills et all 1 19951) or a simplified prescription for heating dur- 
ing the collision was applied to mimic t he presumed pre-mai n 
sequence contraction phase of evolution ( Sandq uist et al.|[T997l) . 
The first realistic evolution calcu lations of collision products 
were performed bv lSills et al.l (1 19971) who took SPH simulations 
of head-on collisions between detailed stellar models and im- 
ported the SPH results into a stellar evolution code. Their mod- 
els demonstrated that none of the previously made simplifica- 
tions are valid: although the collision products are inflated due to 
shock heating, they do not develop substantial convection zones 
during thermal relaxation and do not undergo significant mixing 
during their evolution. 
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This situation changes when the angular momentum of the 
collisions is considered. For collisions that are even slightly off- 
axis, the remnants retain too much angular momentum to re- 
lax into thermal eq uilibrium without reaching break-up veloc- 



of the approximation is the observation that in a stellar model in 
hydrostatic equilibrium the quantity A defined by 



eqi 

ity dLombardi et al.l 19961). T he evolution of such collisions was 
studied bv lSills et alj ( 200 lb . In the absence of a clear mecha- 



nism by which the stars can lose their excess angular momen- 
tum, they artificially removed a large fraction of the angular mo- 
mentum to allow thermal relaxation. Nevertheless, the remnants 
continue to rotate rapidly throughout their main-sequence evolu- 
tion and rotational mixing makes the remnants bluer and brighter 
and significantly extends their lifetimes. The se conclusions w ere 
confirmed by higher-resolution calculations dSills et al.ll2002l) . 

The above studies have focussed on blue stragglers in glob- 
ular clusters, and investigated only a few interesting cases. This 
limitation was imposed by the computation time required and 
numerical difficulties in the evolution calculations. However, the 
importance of stellar collisions for the evolution of star clusters 
calls for a more systematic approach, covering a larger parame- 
ter space and extending to higher masses and younger clusters. 
This is the aim of our current work. We have developed a flexible 
evolution code that is able to evolve stellar collision products un- 
der a wide range of circumstances in a fairly robust manner. As 
our code cannot yet treat rotation and rotational mixing properly, 
for the moment we consider only head-on collisions and ignore 
the effects of rotation. 

As a first step in a systematic study of stellar merger rem- 
nants we have constructed detailed models of the blue strag- 
glers formed by ste llar collisions in the A^-body model of M67 of 
iHurlev et alj (|2005). They evolved a cluster of 36 000 stars from 
zero age to the age of M67 (4 Gyr) taking into account both clus- 
ter dynamics and stellar and binary evolution. In their simulation 
the cluster evolution resulted in 20 blue stragglers at 4 Gyr, eight 
of which had a collisional origin. They formed either as a result 
of dynamical perturbation of a primordial binary, or as a result of 
three-body (binary-single star) or four-body (binary-binary) in- 
teractions. In two of the latter cases, three stars merged in subse- 
quent collisions with the fourth star ending up as a binary com- 
panion to the blue straggler. Hence in total ten collisions were 
involved in the formation of these eight blue stragglers. 

We evolve these collision products with our detailed stel- 
lar evolution code and compare these models with the evolu- 
tion tracks of normal detailed stellar models and fully mixed 
detai led models , as w ell as with the parametric models used 
by IHurlev et alj d2005l) . In particular we investigate the effect 
on the main-sequence lifetime of the merger product (i.e., the 
time during which it will be visible as a blue straggler), its 
position in the Hertzsprung-Russell diagram and the effect on 
the chemical abundances o f the remnant. In a companion pa- 
per dGlebbeek & Polsll2008l) (paper II) we study the influence of 
varying the collision parameters, in particular the masses of the 
two stars and their evolutionary stage. 



2. Tools 

2.1. Modelling the merging process 

To calculate the structure of the collision remnants immediately 
after the collision w e have used the Make Me A Star (MMAS) 
code developed by lLombardi et al.l d2002l) . which produces a 
one-dimensional model that approximates the results of detailed 
smooth particle hydrodynamics (SPH) calculations. The essence 
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(1) 



increases monotonically from the centre to the surface, at least 
as long as radiation pressu re is negligible (See for instance 
lLombardi etaTII 19961 120021 for details) Here, S denotes the en- 
tropy per particle in units of Boltzmann's constant, fi is the mean 
molecular weight of the constituent particles, Bj is the mass of 
particle species i in units of the atomic mass m u , Xj and are 
their abundance fraction by mass and number, respectively. Note 
that A is a function of entropy and composition only. The factor 
a)i denotes the degeneracy of the ground state and is taken as 1 
for nuclei and 2 for electrons. For a classical ideal gas, (Q~|i re- 
duces to 



(2) 



MMAS modifies the initial A profiles of the parent stars to 
correct for shock heating during the collision and estimates the 
amount of mass lost. The amounts of shock heating and mass 
loss are calibrated to detailed SPH results. The remnant profile 
is then constructed by collecting all mass-bins from the parent 
stars, sorting them in order of increasing A and integrating the 
equations of hydrostatic equilibrium and mass conservation. 

The composition profile of the remnant is also determined 
from the A profile. Shock heating leads to partial mixing of ad- 
jacent layers in each of the parent stars. After both A and the 
individual composition profiles have been corrected for shock 
heating, the composition in the remnant is set to the average 
composition of material from the parent stars with that value of 
A. So if for a particular value of A a fraction f\ comes from the 
primary with composition X\ and a fraction fz comes from the 
secondary with composition X2, then the composition X r of the 
corresponding layer in the remnant is 

X t =Xrf l+ X 2 f 2 . (3) 
For more details see lLombardi et ail d2002l) . 

2.2. The stellar evolution code 

To calculate the evolution of the parent stars as well as the fur- 
ther evolution of the collision product we use the stellar evolu- 
tion code (here after STARS) originally developed by Egg leton 
dEggletonlll97ll) and updated by others dPols et al J 1 1 9951) . The 
code uses an adaptive non-Lagrangian mesh that allocates mesh- 
points according to a mesh-spacing function that places more 
meshpoints in regions where a higher resolution is required. This 
means that stars can be evolved with reasonable accuracy with 
as few as 200 meshpoints. 

Our version of STAR S u ses nuclear reaction rates from 
ICaughlan & Fowled dl988l) and|Cayj g hlan et al.ldl98l and opac- 
ities f rom Rog ers & Iglesiasl (1 1992b and lAlexander & Fergusonl 
( 1994). The ass umed heavy-element comp osition is scaled to so- 
lar abundances dAnders & Grevesselll989l) . 

STARS is fully implicit and solves the equations for the 
structure and composition of the star simultaneou sly. Convection 
is tre ated using the mixing-length prescription ( Bo hm- Vitensd 
1958 ) and conv ective mixing is modelled as a diffusion process 
(Egg letorilll972l) . We use a ratio of mixing length to local pres- 
sure scale height II Hp = 2.0. 
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Fig. 1. Entropy (top panel) and hydrogen abundance (lower 
panel) profiles in the merger remnant 2321 (see Table [TJ as a 
function of mass coordinate during the convergence stage. The 
profile at the start of the iteration is indicated by a thin dotted 
line while the target profile is indicated by a thick dotted line. 
The solid line is the final profile after the import is complete. 



Usually the mean molecular weight in stars will be either 
constant in a region or decrease radially outward. In merger rem- 
nants there can be layers in which the molecular weight gradient 
is inverted and a layer of higher mean molecular weight lies on 
top of a layer of lower mean molecular weight. Such a situation 
is unstable and leads to a process know n as thermohaline mixing 
dUlrichlll972HKippenhahn et al.lll980h . The instability arises be- 
cause material with high molecular weight is supported by ther- 
mal buoyancy. When a fluid element exchanges heat with its en- 
vironment it loses buoyancy and begins to sink on the timescale 
for heat exchange, i.e. the local therm al timescale. We mode l 
this as a diffusion process as described in Stancli ffe et al.l d2007). 
The diffusion coefficient is given by the product of the typical 
velocity and size of the fluid elements and the efficiency of mix- 
ing thus depends on their assumed geometry. The efficiency we 
h ave adopted in this work c orresponds to spherical geometry, as 
in Kippen hahn et al] d 19801) . If the fluid blobs are elongated the 
mixing is more efficient. However, since for our adopted choice 
the mixing occurs on the local thermal timescale and is fast com- 
pared to the nuclear timescale, our results are not sensitive to an 
increase in the efficiency of thermohaline mixing. 

STARS uses a si mple model for convective overshooting 
dSchroder et alj[l997l) that allows extra mixing in regions where 
V ra d - V a( j > —6. We have found that this prescription leads to 
spurious mixing in the cores of some of our collision products 
(see 35. IK so we have chosen to disable convective overshooting 
for this work. 



2.3. Constructing starting models for merger remnant 
evolution 

To construct starting models for the collisions we evolved stars 
of the appropriate masses from the zero-age main sequence 
(ZAMS) to the time of collision as listed in Table Q] (see 
These models were then used as input for MMAS. 

The output from MMAS was converted into an input model 
for STARS using implicit calculations starting from a normal 
ZAMS model with the same mass as the collision remnant. This 
model is first evolved until its core hydrogen abundance matches 



that of the remnant. At this point, the entropy and composition 
profiles of the model are adjusted to reproduce the profiles of 
the collision product. Figure Q] shows the profiles at the begin- 
ning and end of the iteration as well as the output profiles from 
MMAS for a generic case. The entropy profile was adjusted by 
adding an artificial energy source e alt to the luminosity equation, 



AL 
Am 



AS 

T— + eart . 



(4) 



Here, e nuc is the net energy generation rate from nuclear re- 
actions, S is entropy per unit mass, T is the temperature and 
-T AS I At is the energy released by gravitational contraction. 
This term vanishes for stars in thermal equilibrium. As men- 
tioned collision remnants are out of thermal equilibrium and the 
effect of e art is to specify TAS /At. A suitable form is given by 



AS 

T— + e art 
At 



AS 
T — , 

T 



(5) 



where AS is the difference in entropy between the current model 
and the target model and t is an artificial timescale on which 
the adjustment is to be made. With this choice an equilibrium is 
reached when the entropy profile in the model matches that of the 
target model, AS = 0. The timescale r is arbitrary in principle 
and can be chosen to change the relative weight of the terms 
in the energy equation, which affects the speed of convergence. 
Choosing r to be of the order of the current timestep was found 
to work well. 

At the same time, the composition profile is changed 
smoothly by setting 



X — > (1 - A)X + AX tMget , 



(6) 



where A is increased monotonically from to 1 in the course 
of the run. Composition changes due to nuclear reactions and 
mixing processes were ignored. We continue this procedure until 
the entropy and composition profiles in the model match those of 
the collision product. The final model is then used as a starting 
model for the evolution of the merger remnant. 

This implicit scheme is very stable and flexible, and deals 
well with small irregularities in the output. Figure Q] shows that 
the resulting entropy and composition profiles agree very well 
with the target profiles. Except for the entropy profile in the core 
the two curves overlap within the thickness of the lines shown. 
We have tested whether small deviations like these affect the 
long-term evolution of the collision product and have found that 
they are unimportant. 

Detailed SPH calculations do not have sufficient resolution 
to resolve the outer parts of the envelope of the collision prod- 
uct and MMAS likewise does not have any real information 
about this region. This means that the structure of the outer enve- 
lope cannot be determined from these models. ISills et al] (1 19971) 
extrapolated the entropy profile and used the condition of hy- 
drostatic equilibrium to reconstruct the outer layers. We have 
found it easier to assume that these layers are in thermal equilib- 
rium and have not tried to impose a particular entropy profile in 
these layers. This is reasonable considering that the local ther- 
mal timescale is short compared to that of the rest of the star. 
Because we do not enforce a large entropy on the o uter layers 
our m odels are somewhat less inflated than those of Sill s et all 
(1997) and consequentl y start at a lower luminosity. As pointed 
out bv lSills et al] ([1997) the long term evolution of the collision 
product is determined by the interior properties and is not sen- 
sitive to the assumptions made for the outer layers. Comparing 
their evolution tracks of collision products to tracks computed 
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Table 1. Initi al conditions for th e collisions and collision product lifetimes. The first column lists the number by which the star is 
refered to in iHurlev et al.l (|2005), followed by the time of collision f, the masses M\ of the primary and M2 of the secondary, the 
total mass M of the remnant according to the detailed models and the remaining main sequence lifetime f ms . For comparison we 
also list the main sequence lifetime T ms of a normal star of mass M, the lifetime T mSj h m of a homogenised star and the remaining 
lifetime of the merger remnant according to the parametric BSE prescription (f mS; BSE)- Time is given in millions of years (Myr) and 
masses are given in solar masses (M Q ). 



ID 


t 


M\ 


M 2 


M 


fms 


Tms 


^"nis,hom 


fms.BSE 


2203 


3480 


1.23 


0.85 


1.93 


539 


1000 


745 


991 


2321 


3960 


1.29 


0.59 


1.74 


599 


1345 


862 


1305 


2565 


3650 


0.95 


0.94 


1.76 


1016 


1297 


1122 


1336 


2973 


3170 


0.89 


0.80 


1.56 


1569 


1850 


1685 


1877 


3121 


3890 


1.09 


0.54 


1.52 


1428 


1993 


1629 


2051 


3289-1 


3610 


0.82 


0.51 


1.25 


3157 


3736 


3736 


3910 


3289-2 


3610 




0.76 


1.86 


1024 


1103 


1058 


1020 


3445 


2770 


0.77 


0.76 


1.41 


2172 


2502 


2361 


2547 


3835-1 


3797 


0.82 


0.60 


1.32 


2657 


3090 


3090 


3191 


3835-2 


3798 




0.31 


1.56 


1765 


1849 


1793 


1776 



with our code for the same masses and metallicity confirms this. 
We find that the contraction timescales are similar and that the 
tracks agree from the main sequence onward. 

In this work we ignore rotation in the collision products by 
assuming that all collisions are head-on. In the more realistic 
case of off-centre collisions the effect of rotation on the evolution 
of collision products can be substantial, as outlined in §Q3 We 
choose to ignore this problem for the moment and we defer a 
discussion of this limitation until jj5.2| 

2.4. The BSE/NBODY4 prescription 

We will compare the outcome of our detailed evolution mod- 
els with the result s obtained in the A^-body calculations of 
IHurlev et al.l d2005l) . These c a lculati on were performed using 
the NBODY4 code lAarsefhl (1 19991) in which binary evolu- 
tion is provided by the Binary Star Evolution (BSE) algorithm 
(IHurlev et al.l l2002). In this algorithm a simple analytic prescrip- 
tion is used to model the outcome of stellar collisions where 
no mass is lost during the collision, i.e. the remnant mass is 
M\ + M2. The merger remnant is replaced by a normal evolved 
main-sequence star with a starting age t' given by 

f ' = i Tn^l Ml — +M 2 — ), (7) 

10 M 1 +M 2 \ r ms ,i T mSi2 / 

so that the remaining lifetime of the collision product is f mSj BSE - 
T ms - f. Eq. (|7]i is based on the implicit assumption that the 
remnant is fully mixed. Here T mSj i, T ms ,i and r ms are the main- 
sequence lifetimes of the two colliding stars and of a normal star 
with mass Mi + M2, respectively. Thes e lifetimes are cal culated 
according to the analytic formulae of IHurlev et al.l d2000l) which 
are based on detailed models that include convective overshoot- 
ing (lPolsetal.1119981) . 

3. Evolution of the merger remnants 

The masses M\ and M2 of the colliding stars and the time of 
collision t a r e take n from the A^-body simulation of M67 by 
IHurlev et all d2005l) and are listed in Table [TJ The collisions 
listed as 3289-2 and 3835-2 are further collisions between a col- 
lision product and another main sequence star. In the case of 
3289-2 this results in a collision product that is more than twice 



as massive as the most massive progenitor star. The most mas- 
sive collision product in our list, however, is 2203. 

We have constructed starting models for the calculation of 
the evolution of the collision products in Table Q] using the 
method described in §|2] All models are calculated for an initial 
hydrogen mass fraction X = 0.70 and a mass fraction of heavy 
elements Z = 0.02. For the double collisions 3289 and 3835 we 
have used the output of the first collision as input for the sec- 
ond collision. All collision products were then allowed to evolve 
until the tip of the giant branch. 

For each evolution model based on the MMAS output we 
have calculated two evolution tracks for comparison. One is a 
star with the same mass as the collision product that is evolved 
from the ZAMS at composition X = 0.70 and Z = 0.02, the other 
is a homogenised version of the collision product. 

Table Q] also lists the total remnant mass M according to 
MMAS, the remaining main-sequence lifetime f ms according to 
our detailed evolution models, the main-sequence lifetime T ms 
of a star with the same ZAMS mass and the lifetime of a ho- 
mogenised version of the collision product, T ms hom- The life- 
times listed for 3289-1 and 3835-1 are the lifetimes these merger 
remnants would have had if they had not been involved in a sec- 
ond collision. Also given in Table [TJ is the remaining lifetime 
according to the BSE prescription, f ms ,BSE ( see §E]f or a compar- 
ison and discussion). 

In Figure|2]we have plotted the evolution tracks of our mod- 
els in the Hertzsprung-Russell diagram. The solid lines are the 
evolution tracks of the detailed collision products while the dot- 
ted lines are the evolution tracks of the normal stars. The ho- 
mogenised models are indicated by a dashed line. We see that in 
general, the homogenised models are hotter and brighter than the 
normal models, while the collision products tend to be brighter 
than the normal stars but less luminous than the homogenised 
models. This is an opacity effect, as will be discussed below. 

3.1. Initial structure and contraction phase 

Initially, all collision products have a large amount of excess 
thermal energy and the star's main energy source is gravitational 
contraction. We refer to this initial phase where the collision 
product is puffed up and out of thermal equilibrium as the con- 
traction phase. 
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Fig. 2. Evolution tracks for the collision products (solid lines) compared to tracks for homogenised models (dashed lines) and a 
main sequence star of the same ZAMS mass (dotted line). Also marked are the positions at 4Gyr for the detailed models (• for 
remnants of a single collision and A for the remnants of two collisions), the homogenised models (0) and the BSE prescription (v). 



There is no significant mixing during the collision. This can 
lead to composition inversions in the remnant: the core is rich 
in helium, on top of which there is a hydrogen-rich layer above 
which there is again a helium-rich layer of material from the core 
of the primary. These composition inversions show up in mod- 



els where the primary star is sufficiently evolved to have burned 
a significant fraction of its central hydrogen to helium while 
the secondary is relatively unevolved, i.e. where the secondary 
is much less massive than the primary. Consequently, they are 
present in models 2203, 2321 and 3121 as well as the two dou- 
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Fig. 3. Composition profiles in the merger remnant 2321, imme- 
diately after the collision (solid line) and on reaching the main 
sequence (dashed line). 

ble mergers 3289 and 3835, although the composition inversions 
are small in this case. As mentioned in section l2~2l a composition 
inversion of this type is secularly unstable and leads to thermo- 
haline mixing. This homogenises part of the central region of the 
remnant. 

As an illustrative example, we will discuss in more detail 
the various evolutionary stages of the collision product 2321. 
This exhibits many of the transient features that occur during the 
evolution of these collision products. 

First we consider the composition profile for collision 2321 
which is shown in Figure [3] During the collision a large por- 
tion of the core of the 1.29M primary has sunk to the centre, 
creating a hydrogen poor core below mass coordinate O.1M . 
Above this region material is slightly mixed with material from 
the O.59M secondary, leading to a hydrogen-rich layer between 
0.2 and O.4M . On top of this there is a helium-rich layer at mass 
coordinate O.6M . Between this layer and the surface of the star 
the composition tends to the primordial composition. 

The helium-rich layer at O.6M is unstable to thermohaline 
mixing, which will mix the helium inward while the star con- 
tracts to the main sequence. At the same time, hydrogen will 
reignite in a shell at O.1M . This burning shell forms as a result 
of a peak in the hydrogen burning rate due to the steep increase 
of the hydrogen abundance in this region. The shell drives a con- 
vection zone that connects to the thermohaline layer and mixes 
the inner O.4M of the star. This mixes helium-rich material into 
this burning shell which has the effect of lowering its efficiency. 
By the time the merger remnant has reached the main sequence 
the burning shell has extinguished. We define the end-point of 
the contraction phase as the moment where central hydrogen 
burning takes over as the main energy source and the star is in 
thermal equilibrium. 

At this point the central hydrogen abundance has increased to 
a mass fraction of 0.46 and the composition profile has changed 
to the dashed line in Figure [3] The collision product has a con- 
vective core on the main sequence that extends to about O.2M . 

The evolution of the core is shown in Figure |H where the 
central temperature T c is plotted against the central density p c . 
The short-dashed lines indicate lines of constant entropy with 
entropy increasing in the direction of increasing temperature and 
decreasing density. The core starts at point a with an entropy that 
is close to that of the core of the former primary. This means 



Fig. 4. Evolution track for the collision product 2321 in a p c - 
T c diagram (thick solid line). Also plotted are the progenitor 
tracks (dashed lines), the track for a star with the same ZAMS 
mass as the collision product (dotted line) and the location of the 
ZAMS (thin solid line). The short-dashed lines represent lines of 
constant entropy for a given composition. The collision product 
starts with an overdense core at the point labeled a and reaches 
the main sequence at point b. 

that the core has an entropy that is too low for a star of its mass 
and it will need to increase its entropy before the star can come 
into thermal equilibrium. This increase in entropy is achieved 
by expansion, which means that the core moves to the left in 
the diagram until it reaches the correct adiabat. At this point the 
core starts to contract again until thermal equilibrium is finally 
reached at point b. At this point the core is slightly hotter than 
the core of a star that was born with the same ZAMS mass. 

3.2. Main sequence evolution 

After the contraction phase the collision products follow more 
or less normal main-sequence tracks, although their luminos- 
ity (Figure 13 and central temperature (Figure HJ are typically 
higher than they would have been for a main sequence star of 
the same ZAMS mass. This is owing to the helium enhancement 
in their envelopes, which increases the mean molecular weight 
p. The luminosity L me i ge r of the collision product scales with the 
luminosity of a normal star ac cording to the homology relation 
dKippenhahn & Weigerlll 19901) 

/ ^merger \ 

^merger mus I | » \°J 

\ ^ms / 

with //merger and yu ms the mass-averaged mean molecular weight 
in the collision product and the main sequence star when these 
have the same effective temperature. The homology relation ^ 
is strictly valid only for homogeneous stars with constant opacity 
but it reproduces the luminosity shift of the merger remnants 
with respect to normal main sequence stars very well. 

The higher effective temperature of the homogenised mod- 
els results from a reduction of the average opacity owing to 
the increased helium content of their envelopes. The opacity af- 
fects the stellar structure by increasing or decreasing the photon 
mean free path. In the weakly mixed collision products, most 
of the helium enhancement is in the compact interior while the 
more extended envelope has the normal ZAMS composition. 
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Fig. 5. Surface [C/Fe] and [N/Fe] abundances as a function of 
the helium core mass (as a measure of age) for collision product 
2321. The carbon abundance (thick solid line) in the collision 
product is depleted when compared to a normal red giant (thick 
dashed line) while the nitrogen abundance (thin lines) is corre- 
spondingly enhanced. 



Fig. 6. Location in the colour-magnitude diagram of the blue 
stragglers at 4Gyr. Detailed models of the single collision rem- 
nants are indicated by • and double collision remnants by A. 
Homogenised models are indicated by O and the BSE prescrip- 
tion is indicated by V. The isoc hrones for t = an d t = 4Gyr are 
constructed from the models by lPols et alJ (Il998h . 



Thus, throughout most of the volume of the star a photon will 
'see' a normal hydrogen rich composition and the opacity is not 
strongly affected. Conversely, in the homogenised models the 
helium enhancement is present in the entire envelope and a pho- 
ton will 'see' a helium-enhanced composition with a lower opac- 
ity, leading to a more compact structure. 

Because the luminosity is enhanced the collision products do 
not lie exactly on the extension of the main sequence, but can lie 
somewhat above it. As a consequence of their higher luminosity, 
the central temperature of the collision products is slightly in- 
creased (see Fig. 2]l and the main-sequence lifetime (i.e. the time 
until core hydrogen exhaustion) is reduced compared to the life- 
time of a normal main sequence star with a similar composition 
in the core. 



3.3. Hertzsprung gap and first giant branch 

We follow the evolution of the collision products through the 
Hertzsprung gap up to the tip of the giant branch. We were un- 
able to calculate models beyond the helium flash. 

The red giant phase is very similar to that of a normal star. 
The first dredge-up occurs when the convective envelope pen- 
etrates into the helium-rich layer. For our example collision 
product 2321 this increases the surface helium abundance to 
T s = 0.334 compared to Y s = 0.296 in a normal red giant of the 
same mass. More of the dredged-up material has undergone pro- 
cessing than in a normal red giant because part of this material 
comes from the core of the original primary. As a result, more 
carbon has been converted to nitrogen by CN cycling, which en- 
hances the nitrogen abundance and depletes the carbon abun- 
dance compared to a normal red giant, both by ~ 0.1 dex in the 
case of collision 2321 as shown in Figure [5] Unfortunately this 
does not exc eed the typical obse rvational error bar of 0.15 dex 
(for instance iGratton et al.ll2000l) . 



3.4. Double collisions 

In the A^-body simulation two blue stragglers are present at 
t = 4 Gyr that resulted from consecutive collisions between 
three stars in a binary-binary interaction. In both cases the sec- 
ond collision happened almost immediately after the first. We 
have not considered the situation where the time between col- 
lisions is long enough for the collision product to have evolved 
before the second collision. We have calculated the outcome of 
the first collision and used that as input for the second collision. 
These double collision remnants evolved in a similar way to sin- 
gle collision products and their main sequence evolution is very 
similar to that of ordinary stars of their mass. This can be un- 
derstood from the fact that their progenitors are all fairly low- 
mass stars that are essentially unevolved at the time of collision. 
For this reason collision product 3289-2 has the highest effective 
temperature (Fig. [2]) and the bluest colour (see Fig|6]in even 
though it is not the most massive collision product. It is likely 
that the structure and evolution of double collisions involving 
more massive and thus more evolved stars, or with more time 
passing between the collisions, will produce remn ants that stand 
out more. The bluest blue straggler in M67 (S977; lMathvsll991l) 
is quite possibly the result of a double collision. 

4. Comparison of different methods 

In Fig. [2] we show as solid symbols the locations in the Hertz- 
sprung-Russell diagram of each of the merger products listed 
in Table Q] after 4 Gyr (th e age of M67 obtained by isochrone 
fitting in lPols et al.lll998l) and compare these to the locations 
of the homogenised models (open symbols) and the positions 
according to the BSE prescription (triangles). In Figure [6] we 
show the corresponding locations for all collision products in a 
colour-magnitude diagram. A comparison between the observed 
blue stragglers and our models is given in paper II. We con- 
verted the theoretical surface parameters (L, T e ff and logg) to 
observable mag nitude and colo ur (My, B - V) using model at- 
mospheres fromlKuruczl (|1992|) with empirical corrections from 
iLejeune et al.l (1 19971) . Also shown is the location of the ZAMS 
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(dashed line) and the t — 4 Gyr isochrone from lPols et al.l (fl998). 
The detailed models are redder than the homogenised models 
and the positions predicted by the BSE recipe, with the largest 
difference occurring for the brightest and most massive blue 
stragglers. 

The homogeneous models are bluer than the detailed models 
for two reasons. First, the helium abundance in the envelopes of 
the homogeneous models is enhanced over that of the detailed 
models, which moves the position of the ZAMS line to the blue. 
Second, if we compare the lifetimes listed in Table Q] we see 
that the remaining lifetime for the homogenised models is longer 
than that of the detailed models, so that the homogenised mod- 
els are closer to their zero-age main sequence. This is particulary 
striking in remnant 2203 (Fig. |2), for which the detailed model 
has almost reached the terminal-age main sequence point. The 
longer lifetime of the fully mixed models is due to the larger hy- 
drogen abundance in their cores compared to the detailed mod- 
els. The largest discrepancy is found for progenitor stars that 
were substantially evolved at the time of collision, in particular 
remnant 2321. The best match between detailed and fully mixed 
models is obtained for the double collision remnants, whose pro- 
genitor stars were not very evolved. In general the approximation 
of homogeneous mixing is worse for collisions involving more 
evolved stars than for collisions involving less evolved stars. 

The BSE recipe gives locations that are brighter and bluer 
mainly because mass loss from the collision is not considered, 
leading to higher remnant masses. As can be seen from Figure|2] 
this effect is most noticeable in remnants 2973 and 3445 and the 
two double collision remnants 3289-2 and 3835-2. Because BSE 
simply uses a normal main-sequence star model for the collision 
remnant, the effect of the higher mass on the luminosity is at 
least partially compensated by neglecting the effect of helium 
enhancement. 

The BSE prescription assumes complete mixing to estimate 
the remaining lifetime of the merger remnant. The BSE lifetimes 
are generally longer than those of both our detailed models (by 
more than a factor 2 in the case of remnant 2321) and the fully 
mixed models, despite the fact that the overestimate of the rem- 
nant mass by BSE would result in an underestimate of the life- 
time by 20-30 per cent. This is, however, counteracted by two 
effects when compared to our fully mixed models. First, the ne- 
glect of helium enhancement leads to an increase of the lifetime 
by up to 30 per cent for collision remnants of evolved progen- 
itor stars. Second, the inclusion of convective overshooting in 
the BSE models gives rise to longer main-sequence lifetimes, by 
25-30 per cent in the mass range considered (see jj5.2| for a dis- 
cussion of this issue). A fairer comparison of the BSE method 
is to compute the remnant lifetime using eq. (O but inserting 
main-sequence lifetimes of non-overshooting models; we then 
find BSE lifetimes that are indeed somewhat smaller than those 
of our fully mixed models. 

The homogeneous models and the BSE prescription can both 
produce blue stragglers that are brighter and bluer than those 
produced by the detailed models. In the case of the BSE models 
this is due to the higher mass of the collision products, while in 
the case of the homogeneous models this is due to the higher 
helium content of their envelopes. This implies that in order to 
form the brightest blue stragglers it may be necessary to collide 
more massive progenitor stars, or that enhanced mixing in the 
collision products (for instance due to rotation) is necessary. 



5. Discussion 

5.1. Convective Overshooting 

There is evidence that for ordinary stars in the mass range we 
consider mixing extends beyond the boundary of the classical 
convective core, which is usually modelled by convective over- 
shooting as discussed in Schrode r et al.l dl997l) . As mentioned 
in section l2~2l we have disabled convective overshooting for the 
evolution of the merger remnants. During the contraction phase 
in some of the remnants, regions develop where V ra( j - V a( j is 
close to zero but remains negative (V, a( j - V a d > -6), so that 
these regions never actually become convective. The overshoot- 
ing prescription used in our code leads to undesired spurious 
mixing in such regions. 

Convective overshooting affects the lifetime and the shape 
of the main-sequence track of a star. The progenitor stars will 
mostly not be overly affected by overshooting (because their 
masses are too small to develop convective cores), whereas the 
collision product is massive enough to experience overshooting. 
The main difference we expect is that for two given parent stars 
the collision product will have a larger reservoir of hydrogen and 
live longer. The inclusion of overshooting in our detailed mod- 
els would likely make the collision remnants of evolved parent 
stars (such as 2203) somewhat bluer, because the enhanced mix- 
ing brings the remnant closer to its ZAMS position. More im- 
portantly, the lifetimes of all our collision products would be in- 
creased, perhaps by about 30 per cent (base d on comparing o ver- 
hooting and non-overshooting models from lPols et al.l l998). 

5.2. Rotation 

We have ignored rotation in the present work by treating all col- 
lisions as head-on. Rotation will modify our results by extra 
mixing and possibly by enhanced mass loss. The latter is par- 
ticularly relevant because off-centre collisions produce remnants 
that have so much angular momentum that their rotation rate will 
approach critical rotation during the contraction phase. How the 
remnant loses angular momentum is still an unsolved problem 
and it is likely that magnetic fields play a role here, possib ly by 
coupling the star to a circumstellar disk dSills et alj l2005h . At 
present, nothing is known about the magnetic field in the colli- 
sion product after the collision. 

Rotational mixing provides a larger hydrogen reservoir for 
nuclear burning and increases the helium content in the enve- 
lope, maki ng the remnants bluer and brighter and extending their 
lifetimes (Sills et al. 2001). By ignoring rotation we implicitly 
assume that an efficient spin-down mechanism operates in the 
collision products such that they avoid both the angular momen- 
tum problem and significant rotational mixing. It is significant 
that the blue stragglers in M67 are slow rotators, with projected 
rotation velocities smaller than ty pical for their spectral type 
dPeterson et al .11 1 9 84t iMath v sll 1 9 9 1 h . Since all formation mecha- 
nisms for blue stragglers are expected to result in rapid rotation, 
this suggests that they can indeed lose angular momentum effi- 
ciently. We will investigate this problem in more detail in future 
work. 



6. Conclusions 

Hydrodynamical simulations of stellar collisions produce rem- 
nants that are out of thermal equilibrium and are not fully 
mixed. We have developed an efficient procedure for importing 
the structural and chemical profiles of such collision products 
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into a fully implicit stellar evolution code. Our evolution code 
is fairly robust and can evolve the collision remnants until the 
tip of the giant branch with a minimum of human intervention. 
We have applied our code to construct detailed models of colli- 
si onal blue stragglers formed in the TY-body simulation of M67 
bv lHurlev et aiT(l2005l) . 



The evolution of collision products depends on the amount 
of mixing during the collision and the thermal relaxation phase. 
Assuming the collision product has been homogeneously mixed 
produces evolution tracks that are too blue while replacing the 
collision product wit h a normal star of the same mass (as done 
in the simulations of iHurlev et al.ll2005l) produces an evolution 
track that is not bright enough. Both approximations overesti- 
mate the lifetime of the collision product. These considerations 
will affect the predicted colour-magnitude diagram distribution 
of collisional blue stragglers from a cluster simulation. 

Our code is suitable for a systematic exploration of the wide 
parameter space of collisions in clusters of different ages. This 
will be the topic of future papers. Eventually we hope to inte- 
grate our code into a full TV-body code to allow for more realistic 
and self-consistent star cluster simulations. 

Acknowledgements. We thank the referee, Alison Sills, for useful comments 
that improved this paper. EG acknowledges support from NWO under grant 
614.000.303 



Appendix A: Modifications to the evolution code 

It is a feature of the STARS code that it solves the stellar struc- 
ture equations simultaneously with the reaction-diffusion equa- 
tions for the different abundances on a moving mesh. The code is 
normally fast and reliable, but we have found that it sometimes 
has difficulty evolving our merger remnants and have come up 
with a scheme that helps it evolve through "difficult" timesteps. 

Consider the set of independent variables H that represents 
a solution to the stellar structure equations at time t. The prob- 
lem is then to find the changes AH such that H + AH repre- 
sents the solution at time t + At. An initial guess for AH can be 
taken from the previous timestep and then improved by iteration 
in a Henyey-like solver. If no solution can be found, a smaller 
timestep can be tried. In some cases, this leads to a runaway sit- 
uation where repeated convergence failures cause the timestep to 
drop until it drops below a threshold value and the code aborts. 

The cause of the convergence failure is that the initial guess 
for the corrections AH is not close enough to the desired correc- 
tions. We have looked for ways to improve the initial guess in 
case of convergence failure. 

Often it is possible to identify the terms in the equations that 
cause difficulty. A common example are the diffusion terms in 
the reaction-diffusion equations, 



dXi 1 , 
at pr z 



(A.l) 



where X t represents the abundance of species i, cr is the sum of 
all diffusion coefficients affecting the composition and Rj is the 
production (or destruction) rate of species i due to nuclear reac- 
tions. In case of convergence problems, it can help te eliminate 
or reduce ("relax") the diffusion coefficient cr. The resulting cor- 
rections AH' are not the final corrections, but they might be a 
better first guess than the values used previously. 

Convective mixing is the most common example where our 
above relaxation scheme is useful, but it is by no means the only 
one. Other examples where we have found it useful in our code 
are the nuclear energy generation rate (which is then relaxed 



from the value at the previous timestep), the mass loss rate and 
advection terms in the luminosity equation. 
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